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FIELD OF THE INVENTION 

The present invention relates to a method intended to obtain reflection travel times 
from an interpretation of seismic data in migrated cylindrical waves, for a given value 
of a parameter defining the slope of these waves, or from the superposition of such data 
5 associated with various substantially parallel acquisition lines, this parameter possibly 
taking successively several values. 

BACKGROUND OF THE INVENTION 

Demigration of seismic data is the inverse of migration, which is a conventional 
seismic data imaging method essentially consisting, knowing the value of a wavefield at 
10 a known depth, for example at the surfece, as well as a distribution model of the 
propagation velocities of the waves in the subsoil, in modelling the propagation of the 
source field and the backpropagation of the reflection data recorded, and in seeking 
phase coherences between these two modelled fields. Migration is particularly usefid 
for interpretation of seismic data acquired on complex structures. 

15 It also allows to access, via demigration, the arrival times of the reflections 

associated with the picked events. This operation is notably carried out for determining, 
by means of kinematic methods such as migration velocity analysis methods, the 

velocity distribution in the subsurface. This determination constitutes a key stage in 
complex structures imaging. An example of such a method is described by : 

20 - Stork, C, 1992; « Reflection Tomography in the Postmigrated Domain »; 
Geophysics, 57, 680-692. 
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A particularly effective migration velocity analysis method, referred to as Smart 
(« Sequential Migration Aided Reflection Tomography ») is notably described by : 

- Ehinger, A. and Lailly, P., 1995, Velocity Model Determination by the SMART 
Method, Part 1 : Theory : 65* Ann. Intemat. Mtg., Soc. Expl. Geophys., Expanded 

5 Abstracts, pp. 739-742, 

Reflection tomography is used, in such methods, for determining the velocity 
model ; the reflection travel times (reflection tomography data) are obtained, in the 
complex parts, by demigrating the picked events upon interpretation of the migrated 
data. This method has been successfully implemented for 3D imaging of complex 
10 structures by using 2D prestack migration algorithms (and the associated demigration) 
on a set of seismic profiles extracted from a 3D acquisition. However, such an approach 
is efficient only in a very specific context : 

- on the geologic plane, the structure must vary relatively slowly in a direction 
(referred to as « strike » in the Anglo-Saxon literature) ; 

15 - the seismic profiles must be acquired in the orthogonal direction (therefore the 
complexity direction) referred to as « dip » in the Anglo-Saxon literature. 

On the other hand, this approach is difficult to use with complex structures without 
a privileged direction : such structures require prestack 3D migration techniques, which 
involve long calculation times, especially if imaging is carried out on a volume. 

20 Patent FR-2,784,195 filed by the applicant notably describes another seismic data 

migration method known to the man skilled in the art as Alchemig. This method allows 
to perform a 3D prestack migration of seismic events and imaging of the volumes of an 
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underground zone, with an attractive calculation time, from a series of a number Ns of 
reflection survey cycles. 

It comprises successive emission of elementary wavefields, each defined by the 
association of a seismic signal W(t) and of a determined emission point from a series of 
5 emission points 5, with 1 < i < N^, reception, by seismic receivers placed at positions 
-R/ , of the seismic signals reflected by the zone in response to each one of these 
wavefields, and recording of the various signals received by each seismic receiver in 
form of time-dependent seismic traces dj(f) . For a given velocity model, 

a) a slowness vector p (homogeneous to the inverse of a velocity) whose two 
10 con^nents py and px can take a previously defined series of values is defined, 

b) a time lag ftmction fg(p,i) is defined for a given slowness vector p and for a given 
emission point <S, , 

c) time lag ftmction t^ipj) is applied to each elementary wavefield associated with 
emission point S, and a conq)osite wavefield is formed at the surfece by spatioten^ral 

15 superposition of the various elementary wavefields to which such a time lag is applied, 

d) a time lag t^(p,i) is applied to each seismic trace df(t) marked by pair (ij) and a 
conqjosite trace field is formed at the surfece by spatiotemporal superposition of the 
various seismic traces to which such a time lag is applied, 

e) a migration of the conposite trace field is carried out by using as the wavefield the 
20 conqjosite wavefield, which is done by modelling the propagation of the composite 

wavefield and the baclqpropagation of the composite trace field, and by suitably 
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combining the two composite fields thus modelled at ^ly point of the zone to be 
imaged, 

f) stages c) to e) are repeated for all the values taken by components Py and px of vector 

p , and 

5 g) forany fixed value of the second con5)onentp^ of said vector p , the results of these 
various combinations are summed so as to obtain a n%rated image associated with this 
fixed value of p„ thus performing a prestack migration. 

It is assumed that, even if the technique proves robust in case of violation of these 
hypotheses, the seismic data is constant-azimuth data (the azimuth being defined as the 
10 direction of the source-pickup bipoints) (hypothesis 1) and that acquisition is carried out 
by displacing the source along lines parallel to the direction defined by the azimuth 
(hypothesis 2). 

If this is not the case, it can be remedied by first carrying out a preprocessing 
operation, for example according to the technique referred to as AMO (« Azunuth 
15 Moveout »), well-known to the man skilled in the art, to adjust the data obtained from 
standard marine acquisitions. 

The co-ordimtes sj^em selected here is such that the x-axis is parallel to the 
direction of the acquisition lines and the y-axis thus represents the position of an 
acquisMon line. We define hereafter by « cylindrical wave data » the measurement, 
20 along an acquisition line, of the seismic response of the subsurfece to a cylindrical wave 
(wave generated by a source line and whose phase varies linearly with the abscissa 
along the acquisition line) whose axis comcides with the acquisition line considered. 
Considering this definition, the Alchemig technique allows to obtain the superposition. 
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along the various acquisition lines, of « migrated cylindrical wave data », these 
cylindrical waves being associated with an acquisition line and with a predetermined 
value of parameter Px defining the slope of the cylindrical wave. The procedure consists 
in applying to the data linear phase shifts parametered by vector p whose component 
5 along the profiles axis is px, and in calculating such a superposition by applying to the 
phase shifted data a plane-wave migration software. Postmigration stacking is obtained 
by summing the results obtained for the different values taken by p^. We suppose that, 
by suitable selection of the origio of the mark or of the phase shift fimction, the point of 
the surface where the phase shifts are zero has (x = 0, y = 0) as the co-ordinates. 

10 It can be noted that the technique described above is one way among many to obtain 

superposition of migrated cylindrical wave data associated with parallel acquisition 
lines and with a given value of parameter p^. Another procedure can for example 
consist, as described by Claerbout, J.F., 1971, « Towards a Unified Theory of Reflector 
Mapping » ; in Geophysics, 36, No.3, 467-481, in calculating the result by carrying out, 

15 for each acquisition line : 

- propagation of the cylindrical wave and backpropagation of the cylindrical wave 
data, and 

- calculation, at any point of the space, of the crosscorreiation of the source field 
whose propagation has been modelled and of the backprop^ated data field, which 

20 allows to obtain, for the acquisition line considered, the migrated cylindrical wave data ; 

and by stacking (superposing) the migrated cylindrical wave data obtained for the 
various acquisition lines. 
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If we have, for several values of parameter p„ the stq)erposition of migrated 
cylindrical wave data associated with various acquisition lines, it is possible, as 
explained in the aforementioned patent, to control the quality of the velocity model used 
for the migration and also to update the velocity model by exploiting the deformations 
5 of the events when going from one cylindrical wave parameter to the next (migration 
velocity analysis). 

For updating, according to the aforementioned Smart method (or others), the events 
have to be picked, for various values of parameter p„ in the superposition of the 
migrated cylindrical wave data, a superpoation corre^nding to the various acquisition 
10 lines, then the picked events have to be demigrated and the model is updated by 
reflection tomography. Even so one needs to know how to demigrate the interpreted 
migrated events when interpretation has been performed on migrated cylindrical wave 
data or their superposition obtained for various acquisition lines. 

SUMMARY OF THE INVENTION 

15 This is one of tl:^ con^lex problems that can be solved by Hbe method according to 

the invention. The approach selected furthermore leads to methods that will be 
described, and which are closely related to the concept of cylindrical wave. 

The method according to the invention allows to determine the reflection travel 
times of seismic events picked on 3D records of seismic data corresponding to signals 
20 picked up by seismic receivers distributed along an acquisition line (respectively along 
several acquisition lines) in response to the emission, in the subsoil, of waves from one 
or more source points, this data having been previously converted into migrated 
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cylindrical wave data and intapreted (respectively in a superposition of migrated 
cylindrical wave data and interpreted). 

In cases where interpretation has been carried out on the migrated data associated 
with a single cylindrical wave, and therefore a single acquisition line, the method 
5 comprises the following stages : 

a) defining a slowness vector {p) whose component (p^) ia a direction parallel to the 
acquisition line defines the slope of the cylindrical wave ; 

b) for a seismic receiver situated at abscissa (Xr) on the acquisition line, seeking 
abscissa (^) of the source point on the acquisition line such that a ray starting from the 

10 seismic receiver and reflecting on a picked event emerges at said source point, with a 
slowness vector (p) whose con^onent in the direction of the acquisition line is (-px) ; 

c) determining a travel time (t^ (Xr)) by adding to the value of the travel time along the 
ray a time equal to the product of the slope of the cylindrical wave by the abscissa of the 
source point ; 

15 d) repeating stages b) and c) for all the positions of the receivers for which a 
demigration result is wanted ; and 

e) repeating stages a) to d) for all the values taken by parameter (p^) for which an event 
has been picked in the corresponding migrated data. 

In cases where ioterpretation has been carried out on the superposition of 
20 cylindrical wave data, a superposition performed from several acquisition lines, the 
method comprises the following stages : 



a) defining a slowness vector Cp) whose conqwnent (p^ in a direction parallel to the 
acquisition liaes defines the slope of the cylindrical waves associated with the various 

acquisition lines ; 

b) for a given seismic receiver situated at abscissa (Xr) on an acquisition line, seeking 
5 abscissa (4) of a source point on the acquisition line such that a ray starting from the 

seismic receiver and reflecting on a picked event emerges at said source point, with a 
slowness vector (p) whose conponent in the direction of the acquisition line is (-p^ ; 

c) determining a travel time (f (Xr)) by adding to the value of the travel time along said 
ray a time equal to the product of the slope of the cylindrical wave by the abscissa of the 

10 source point ; 

d) repeating stages b) and c) for all the positions of the receivers for which a 
demigration result is wanted ; 

e) repeating stages a) to d) for all the acquisition lines for widch a demigration result is 
wanted, 

15 f) repeating stages a) to e) for all the values taken by parameter (p^ for which an event 
has been picked in the corresponding migrated data. 

The prior processing operation leading to migrated cylindrical wave data can 
possibly involve, if necessary, AMO type preprocessing as mentioned above or a spatial 
resampling of the seismic data acquired in the field. 

20 According to an implementation mode, the travel times associated with the 

cylindrical wave events are used to implement a cylindrical-wave reflection tomography 
technique eaxd to calculate a velocity distribution in the environment. 
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According to another implementation mode, the travel times associated with the 

cylindrical wave events are converted into travel times associated with the events by 
shotpoints, by exploiting the correspondences between the cylindrical wave data and the 
shotpoint data 

5 According to an implementation mode, when the seismic data acquisition device is 

organized by acquisition lines, le. when a series of source points and a series of pickups 
are arranged in one line, the synthesized shot travel times are converted to travel times 
associated with the source-receiver pairs (of the seismic data acquisition device), this 
conversion being obtained by extrapolation in the direction of the acquisition lines of 
10 the synthesized travel time for the source-receiver pair that is the closest to the source- 
receiver pair for which information on the travel time is wanted. 

According to an implementation mode, the seismic data being only constant- 
azimuth data, the synthesized shot travel times are converted into travel times associated 
with the source-receiver pairs, this conversion being performed by means of a double 
15 extrapolation, one in the direction defined by the azimuth, ibe second in the direction 
orthogonal to the previous direction. 

The method can be used notably to determine the velocity model by applying for 
exan^le a reflection tomography to the travel times obtained. 

BRIEF DESCRIPTION OF THE FIGURES 

20 Other features and advantages of the method accordii^ to the invention will be 

clear from reading the description hereafter of non-limitative examples, with reference 
to the accon:panying drawings wherein : 
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- Figure 1 describes the various options allowing implementation of the method 
according to the invention, and 

- Figure 2 compares, on a constant-offset record section, the seismic data with the 
arrival times obtained by implementing an option of the method. 

5 DETAILED DESCRIPTION 

Demigration of events picked on migrated cylindrical wave data requires a specific 
technique. We will successively distinguish : 

- demigration of an event picked on migrated cylindrical wave data, 

- demigration of the superposition (stack) of migrated cyhndrical wave data, 
10 associated with various acquisition lines, 

- obtaining kinematic data (reflection travel time) associated with source-pickup 
pairs, and 

- obtaining kinematic data associated with source-pickup pairs existing in the original 
acquisition device (Le. for the source-pickup pairs used for migration). 

15 1) Demigratfon of an event picked on migrated cylindrical wave data 

The result of such a demigration is the time of arrival, on the records associated 
with the cylindrical wave considered (characterized by the co-ordinate along the y-axis 
of the acquisition line and by the value of parameter px definii^ the slope of the 
cylindrical wave), of the reflection corresponding to the picked event. This arrival time 

20 t* depends on the abscissa Xr of the receiver considered. Demigration is naturally 
conducted by seeking, for a given receiver position, the abscissa E, of the source point or 
of the sotirce points S on the line considered, such that the ray starting from the receiver 
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and reflecting on the picked event emerges at S with a slowness vectcff p whose 
conponent along x is -px, the velocity model beii^ identical to the model used for 
migration. 

Thus worded, the problem appears to be the solution of a non-linear equation 
(e5q>ressing the condition on the slowness vector) in a single unknown (abscissa % of the 
source). This equation can be solved numerically by means of any method suited for 
solution of non-linear equations, such as Newton's method for exan^le. 

The travel time resultii^ from the demigration t* (Xr) will be obtained by taking the 
value of the travel time along all the rays obtained for source S, a time to which quantity 
px^ is added. It can be noted tiiat th^e is in fact an infinity of such rays, but these rays 
give, for a given arrival, the same travel time. This operation is repeated for all the 
receiver abscissas for which the demigration result is sought. The arrival times thus 
obtained constitute intrinsic inform^ition connected with the synthesized cylindrical 
wave data. They do iwt depend on the velocity model used for tiie migration (which 
justifies the feet that they can be used for determination of the velocity model), because 
the demigration operation described above is, as regards kinematics, the inverse of the 
migration of cylindrical wave data. Component Pr along x of slowness vector p of the 
ray at the emergence of the seismic pickup also constitutes a quantity that is intrinsically 
linked with the data and vviiich therefore does not depend on the velocity model used for 
the migration ; it gives the slope, on the cylindrical wave data, of the event recorded by 
the pickup. 
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2) Demigration of the superposition of migrated cylindrical wave data, associated 
with various acquisition lines 

This operation is simply carried out by repeating the operation defined in the 
previous paragraph for all the acquisition lines where a demigration result is desired and 
5 it delivers the same intrinsic quantities, i.e. the time of arrival of the event, on the 
cylindrical wave data associated with each line, as a fimction of the abscissa of the 
pickup and of the slope of this event on the x-axis, i.e. component pR along x of the 
slowness vector when the ray emerges at the pickup. Furthermore, if we denote the 
components on the y-axis of the slowness vector for rays emerging at the surface 
10 respectively at the pickup and at the soxirce by and , quantity + pi also 
appears as an intrinsic quantity ; it is interpreted as the slope, on the y-axis (therefore 
from one line to the next), of the event seen for the seismic pickup abscissa considered, 
on cylindrical wave data associated with parameter p^. 

3) Obtaining kinematic data (reflection travel time) associated with source-pickup 
15 pairs 

If a reflection tomography software allowing to process travel times associated with 
cylindrical waves is available, the data obtained in the previous paragraph is sufficient 
for updating the velocity model and for migration velocity analysis according to the 
Smart method. 

20 If such inversion softwares are not available, and if conventional softwares 

processing kinematic data organized by source-receiver pairs are used, the stage 
described hereafter will be carried out. 
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This stage consists in exploiting correspondences between the cylindrical wave data 
and the shotpoint data, Le. all the records obtained with the various receivers for the 
same seismic source assumed to be a point source. 

Such a correspondence can be exploited only if cylindrical wave data associated 
5 with several values of parameter p^ defining the slope of the cylindrical waves is 
available. It can be noted that migration has no role here : the point is simply to 
transform data. The kinematic data (reflectk>n travel time) associated with source- 
pickup pairs is obtained by proceeding acquisition line by acqidsition line (the source- 
pickup pairs are assumed to be situated on these lines). 

10 For a given line, we have, by applying the techniques described in the previous two 

paragraphs, function r (Xr) for various values of parameter p^. In other words, we have 
the fimction, that we shall call, after a notation change, r (p^, Xr), or more exactly a 
sampling of this function, the values taken by this fimction being obtair^d for specific 
values of the (p« Xr) pair. A fimction Xs (p« Xr) is defined by the relation : 



15 a quantity that we can evaluate nimierically (for example by using approximations 
based on finite differences) for specific values of the (p„ Xr) pair. Formula (2) hereafter 
gives the time denoted by ()^,Xr) of the event at the pickup of abscissa Xr in the 
shotpoint data for a shot abscissa Xs = ^ (p» Xr). 




(1) 



(2) 
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The times thus obtained are intrinsic quantities : they are only related to the seismic 
data and they do not depend on the way migration has been implemented. In conclusion, 
it is sufficient to ^ply formula (2) on the various lines for all the (p„ pairs for which 

data t* (p„ Xr) and ^— (p^, Xg) is available to know the arrival times of the interpreted 

event (or, in practice, of each interpreted event) for a multitude of source-pickup pairs 
covering the surfece and thus allow precise determination, via a reflection travel time 
tomography, of the velocity model. 

4> Obtaining k inematic data associated with source-pickup pairs existing in the 
original acquisition device 

Even if tte data obtaii^ by means of the technique described in the previous 
section is sufficient to update the velocity model, it may be useful to compare the times 
obtained with the non-migrated original seismic data. Such a conparison is in feet 
interesting for validating the picking (and notably the interpretation) performed on the 
migrated data and for checking the reliability of the demigration which, as mentioned 
above, involves numerical procedures based on approximations and ray tracing, a very 
sensitive tool for numerical representation of the reflectors. 

However, direct comparison is difficult since there is no reason for the source- 
pickup pairs to be those corresponding to the data used for the migration. In feet, the 
travel times provided by the method described in the previous section are obtained for 
source-pickup pairs distributed as a series of clusters distributed in 2D spaces (each one 
of these spaces corresponding to an acquisition line), whereas the original data can be 
distributed in form of a cluster distributed in a 3D space. However, even in cases where 
migration has been carried out fix)m data organized by 3D lines, the problem remains 
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the same : there is no reason for the source-pickup pair defined by the pair of abscissa 
(Xs (p„ Xr ), Xr), obtained as explamed in the previous section, to appear in the original 
acquisition device. 

Interpolation or extrapolation techniques are necessary to obtain kinematic data 
5 associated with source-pickup pairs existing in the origioal acquisition device from the 
results obtained according to the method described in the previous section. One 
technique consists in defining a proximity criterion in a 2D or 3D space and, for a given 
source-receiver pair of the original acquisition device, in assigning thereto the travel 
time of the source-pickup pair that is the closest among those for which travel time 
10 information has been obtained by means of the method described in section 3. Such an 
approach is viable if the information extrapolation distance is not too great. In the 
opposite case, more elaborate techniques have to be implemented, such as those 
described hereafter, according to whether resan^led data or constant-azimuth data is 
used. 

15 4-1 Implementatbn on data organized by acquisition lines 

In this implementation type, we also suppose that the seismic data is constant- 
azimuth data and that acquisition is carried out by displacing the source along lines 
parallel to the direction defined by the azimuth, for exanq)le after spatial resampling of 
the data and/or AMO type preprocessing, as mentioned above. The source-pickup pairs 
20 involved in the implementation of the migration method are distributed as clusters 
distributed in 2D spaces, each one of these spaces corresponding to an acquisition line. 

In this case, the data is organized by acquisition line and we want to obtain the time 
of arrival of the event reflected for the source-pickup pair of the original acquisition 
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device (disp), defined by the acquisition line to vMch it belongs and by abscissas 
xf^ and xf^ . As explained above, we have extracted from all the arrival times 
calculated for the line considered the arrival time associated with the source-pickup pair 
defined by pair (xf'',xf') that is the closest to the time associated with pair 
5 (xf'^,xf^) among the cluster obtained in the 2D space described by pairs (Xs, Xr). We 
then use formula (3) given hereafter to obtain the arrival time t^{xf^^,x^^) from the 
quantities obtained during the previous processing stages : 

f^xf ^xf O = t\xr'',xf^) - pAxt" - xf^) + p^dxp' - xf 0 (3) 

quantity Pr having been defined in sections 1, 2 above. It can be noted that this formula 
exploits the following important result : according to hypotheses 1 and 2 stated in 
10 section I, for a given pickup and therefore for a given acquisition line, the sk)pe of an 
event in the cylindrical wave data associated with parameter p^ is identical to the slope 
of the event in the shotpoint collection for the shot related by relation (1) to pair (p» Xr) 
of abscissa (x^ of the pickup. The aforementioned extrapolation technique is of course 
not limitative. An interpolation technique can be used just as well. 

15 4-2 Implementation of the method on constant-azimuth data 

In this implementation type, we assume that the seismic data is constant-azimuth 
data. The source-pickup pairs involved in the implementation of the migration are 
distributed as clusters distributed in a 3D space, and these pairs can be described by 
means of three parameters : the respective abscissas of the source and of the pickup, and 
20 their common ordinate. The problem consists here in obtaining the arrival times of the 
event reflected for the source-pickup pairs of the original acquisition device (therefore 
distributed as clusters in a 3D space) from the calculated times, by using the method 
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presented in section 3. This metliod is tben implemented for a SCTies of fictitious 
acquisition lines of ordinate yi, i representing the number of the line considered. These 
lines have to be sufficiently close to one another to allow extrapolation, on the source- 
pickup pairs of the original acquisition device, of the travel times calculated according 
to the method described in section 3. We thus want to calculate the travel time of the 
reflected event for a source-pickup pair of the original device, this pair having 
(x^.xp') as the abscissa and y*^ as the common ordinate. 

We therefore extract, from all the arrival times calculated for the series of lines 
mentioned above, the time associated with the source-pickup pair defined by pair 
x^%x^) and by ordinate y; that is the closest, according to the metric defined in the 
3D space, to the time associated with the source-pickup pair considered of the original 
device. We then use formula (4) hereafter to obtain the arrival time t\xf^,x^^'',y'^'^) 
from the quantities obtained during the previous processing stages : 

nxr^xp'y^) = t\xf\xf^,y,)-p^(xt^-xt) + p,(xt'-xr^) + (py + py^^^^^^ 

Here again, the extrapolation technique is given by way of non limitative example. 
An interpolation technique can be used just as welL 

Figure 2 con5>ares, on a constant-offset sectk>n (1500 m here) obtained after 3D 
seismic prospecting operations in the North Sea, the seismic data with the arrival times 
obtained by implementing the technique described above. The reflection can be clearly 
seen at the base of the salt structure. 



